#include "stdafx.h" 
#include "sequence.h" 
#include "AlignSequence.h"
#include "timer.h" 

using namespace Concurrency; 


bool AlignSequence::CreateLibrary (Sequence & Library) //'::' tells that its two different methods, function that belongs to the class
{
	const int hashLength = 10; 
	int hashTableLength = 1 << (2*10); //this is binary equivalent to 2^10 
	int* hashTable = (int*)malloc (hashTableLength * sizeof(int)); //need enought room for hash table length intergers, getting memory (malloc) 
	memset (hashTable, 0, hashTableLength * sizeof(int));
	char* pLib = Library.sequenceData;
	int collsionCount =0;
	for (int seqLoc = 0; seqLoc < Library.sequenceLength; seqLoc++)
	{
		int HashVal =0;
		for (int j = 0; j < hashLength; j++)
		{
			char Nucleotide = pLib [seqLoc+j];
			int NucValue = 0; 
			switch (Nucleotide)
			{
			case 'A': 
				{
					NucValue = 0;
				}
				break;
			case 'G': 
				{
					NucValue = 1;
				}
				break;
			case 'T': 
				{
					NucValue = 2;
				}
				break;
			case 'C': 
				{
					NucValue = 3;
				}
				break;
			default:
				break;
			}

			HashVal = HashVal << 2 ; //The '<<' operator will shift bits to the left when the char on the left is a number and not a 'cout' 
			HashVal += NucValue;
		}//end of Jloop

		//nucValue is the index (hashValue) of the sequence. we store the location of it (i) at this index location
		if (hashTable [HashVal]!=0)//hash Collision! 
		{
			collsionCount ++;
		}
		hashTable [HashVal] = seqLoc; 
	}
	return (true); 

}



AlignSequence::AlignSequence(void)
{
	secondsToAlign = 0.0; 
}


AlignSequence::~AlignSequence(void)
{
}

double AlignSequence::getTimetoAlign ()
{
	return secondsToAlign; 
}


//int AlignSequence::AlignAMP (Sequence & Probe, Sequence & Library)
//{
//	return 0;
//}

int AlignSequence::AlignAMP (Sequence & Probe, Sequence & Library,int nucleotidesPerHash )
{
	Timer timer;
	timer.Start ();

		//create hash tables for probe and library
	Probe.NucleotiesInHash = nucleotidesPerHash;
	Library.NucleotiesInHash = nucleotidesPerHash;



	//Concurrency::extent<1> libraryExtent ((Library.sequenceLength/4)-3);
	Concurrency::extent<1> libraryExtent (Library.sequenceLength-12);// need to round number of nucleotides in hash to a multiple of 4 
	Concurrency::array_view<const unsigned int, 1> libraryArrayView (Library.sequenceLength/4,(unsigned int*) Library.sequenceData); // the array view handles copying data before the cPU memory and the GPU memory 

	int HashTableLength = Library.getHashTableSize();
	Concurrency::extent<1> libraryHashExtent (HashTableLength);
	Library.allocateHashTable(); 

	//Compare performance for int and unsigned int and const unsigned int
	Concurrency::array_view<unsigned int, 1> libraryHashArrayView (HashTableLength,(unsigned int*) Library.hashTable); // the array view handles copying data before the cPU memory and the GPU memory 


	//Concurrency::extent<1> probeExtent ((Probe.sequenceLength/4)-3);
	//Concurrency::array_view<const unsigned int, 1> probeArrayView (Probe.sequenceLength/4,(unsigned int*) Probe.sequenceData); // the array view handles copying data before the cPU memory and the GPU memory 

	int NucleotidesPerHash = Library.NucleotiesInHash;
	libraryHashArrayView.discard_data();


	parallel_for_each(libraryExtent, [=] (index<1> idx) restrict (amp) {
	int libraryOffset = idx [0];  // current threads index from libraryExtent, zero to library size 

	///create hash for libraryOffset
	
	  // Calculate the hash for first pattern_length characters of the record (exclude null terminating characters)
// Copyright (c) Microsoft Corp. 
//
// File: ampstrstr.cpp
//
// Licensed under the Apache License, Version 2.0 (the "License"); you may not use this 
// file except in compliance with the License. You may obtain a copy of the License at 
// http://www.apache.org/licenses/LICENSE-2.0  
//  
// THIS CODE IS PROVIDED *AS IS* BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, 
// EITHER EXPRESS OR IMPLIED, INCLUDING WITHOUT LIMITATION ANY IMPLIED WARRANTIES OR 
// CONDITIONS OF TITLE, FITNESS FOR A PARTICULAR PURPOSE, MERCHANTABLITY OR NON-INFRINGEMENT. 
//  
		//This hash function assumes that the characters are encoded four each into an int array
	//  Offest is still in chars
	//  int1   int2  int3
	//[ ACGT] [CGTA] [GTAC]

	//First, we need to figure out how walk bytes in the stream starting at offset
	// Divide by 4 will get the correct integer -- so for example, if we start at byte 2, 2/4 = 0 in integer division
	// we use the modulo (remainder to figure out how much to shift
	

	
		int HashVal =0;
		for (int j = 0; j < NucleotidesPerHash ; j++) //Cannot use Library.NucleotiesInHash because Library object uses Chars
		{
		// Divide by 4 will get the correct integer -- so for example, if we start at byte 2, 2/4 = 0 in integer division
		unsigned int cval =  libraryArrayView ((libraryOffset + j)/4);
		// we use the modulo (remainder to figure out how much to shift
		int remainder = (libraryOffset + j) % 4;

		//Integers are stored in backwards -- least significant byte first  (ACGT in memory corresponds to TGCA as a 32 bit integer)
		// if the remainder is 0, we are at 0, 4, 8, etc    ACGT so we want A
		//			which encodes least significant byte  so we need to shift to right 3*8 bits
		// if the remainder is 1, we are at 1, 5, 9, etc    ACGT so we want C,
		//			so we need to shift to right 1*8 bits
		// if the remainder is 2, we are at 2, 6, 10, etc   ACGT so we want G, so we need to shift to right 2*8 bits
		// if the remainder is 3, we are at 3, 7, 11, etc   ACGT so we want T, so we need to shift to right 3*8 bits
		cval = cval >>  remainder *8;
	
		//now we need to mask off the upper three bytes by ANDing with 0x000000FF;
		cval =cval & 0x000000FF;

	
			int NucValue = 0; 
			switch (cval)
			{
			case 'A': 
				{
					NucValue = 0;
				}
				break;
			case 'G': 
				{
					NucValue = 1;
				}
				break;
			case 'T': 
				{
					NucValue = 2;
				}
				break;
			case 'C': 
				{
					NucValue = 3;
				}
				break;
			default:
				break;
			}

			HashVal = HashVal << 2 ; //The '<<' operator will shift bits to the left when the char on the left is a number and not a 'cout' 
			HashVal += NucValue;
		}//end of Jloop

		

       
		libraryHashArrayView (HashVal) = libraryOffset; 
	///end create hash 


	
	});//end of parallelforeach 
	libraryHashArrayView.synchronize (); 

	int* LibraryList = (int*)malloc (Probe.sequenceLength * sizeof(int));	
	memset (LibraryList, 0, Probe.sequenceLength  * sizeof(int));
	double averageOffset = 0.0;
	int probeHitCount = 0; 

	//create hash tables for probe and library
	Probe.createHashTable ();
	//Library.createHashTable (); 
	for (int hashValue = 0; hashValue < Probe.getHashTableSize(); hashValue++)
	{
		int ProbeLocation = Probe.hashTable [hashValue]; 
		if (ProbeLocation == 0)
		{
			continue; 
		}
		//does library have this value?
		int LibraryLocation = Library.hashTable [hashValue]; 
		if (LibraryLocation == 0)
		{
			continue; 
		}
		//we found a match 
		int AlignmentOffSet = (LibraryLocation - ProbeLocation); //how far is the probe found into the library 
		LibraryList [ProbeLocation] = AlignmentOffSet; 
		//averageOffset = averageOffset + AlignmentOffSet; 
		/*probeHitCount ++; 
		cout << " \n Alighnment Off Set = " << AlignmentOffSet << " Average = " << averageOffset/probeHitCount << endl; */

	}
	std::sort(&LibraryList[0], &LibraryList[Probe.sequenceLength-Probe.NucleotiesInHash-1]);

	int commonValue = 0; 
	int commonValueCount = 0;
	int maxCommonValueCount =0;
	int maxCommonValue =0;
	int zeroCounter = 0; 

	for (int i = 0; i < Probe.sequenceLength - Probe.NucleotiesInHash; i++)
	{
		int value = LibraryList[i] ;
		if (value == 0)
		{
			zeroCounter++; 
			continue; 
		}
		if (commonValue == value)
		{
			commonValueCount++;
			if (commonValueCount > maxCommonValueCount)
			{
				maxCommonValueCount = commonValueCount;
				maxCommonValue = value; 
			}
		} 
		else
		{
			if (commonValueCount > 2)
			{
				//cout << "Found " << commonValue << " a total of " << commonValueCount <<  " times\n";
			}
			commonValue = value; 
			commonValueCount = 0; 
		}
	}

	//maxCommonValue is the result to return
	free (LibraryList);
	timer.Stop (); 
	secondsToAlign = timer.Elapsed (); 
	cout << "total library misses =" << zeroCounter << endl;

	return maxCommonValue;
}



int AlignSequence::InstrumentedAlign (Sequence & Probe, Sequence & Library, int nucleotidesPerHash)
{
	Timer timer;
	timer.Start ();
	
	//create hash tables for probe and library
	Probe.NucleotiesInHash = nucleotidesPerHash;
	Library.NucleotiesInHash = nucleotidesPerHash;
	Library.InstrumentedcreateHashTable (); 
	int returnVal =  Probe.ResolveAlignment(Library);
	
	timer.Stop (); 
	secondsToAlign = timer.Elapsed (); 
	return returnVal;

}//end of 


int AlignSequence::FastAlign (Sequence & Probe, Sequence & Library, int nucleotidesPerHash)
{
	Timer timer;
	timer.Start ();
	
	//create hash tables for probe and library
	Probe.NucleotiesInHash = nucleotidesPerHash;
	Library.NucleotiesInHash = nucleotidesPerHash;
	Library.createHashTable (); 
	int returnVal =  Probe.ResolveAlignment(Library);
	
	timer.Stop (); 
	secondsToAlign = timer.Elapsed (); 
	return returnVal;

}//end of 





int AlignSequence::Align (Sequence & Probe, Sequence & Library, int nucleotidesPerHash)
{
	Timer timer;
	timer.Start ();
	
	//create hash tables for probe and library
	Probe.NucleotiesInHash = nucleotidesPerHash;
	Library.NucleotiesInHash = nucleotidesPerHash;


	int* LibraryList = (int*)malloc (Probe.sequenceLength * sizeof(int));	
	memset (LibraryList, 0, Probe.sequenceLength  * sizeof(int));
	double averageOffset = 0.0;
	int probeHitCount = 0; 

	
	Probe.createHashTable ();
	Library.createHashTable (); 
	for (int hashValue = 0; hashValue < Probe.getHashTableSize(); hashValue++)
	{
		int ProbeLocation = Probe.hashTable [hashValue]; 
		if (ProbeLocation == 0)
		{
			continue; 
		}
		//does library have this value?
		int LibraryLocation = Library.hashTable [hashValue]; 
		if (LibraryLocation == 0)
		{
			continue; 
		}
		//we found a match 
		int AlignmentOffSet = (LibraryLocation - ProbeLocation); //how far is the probe found into the library 
		LibraryList [ProbeLocation] = AlignmentOffSet; 
		//averageOffset = averageOffset + AlignmentOffSet; 
		/*probeHitCount ++; 
		cout << " \n Alighnment Off Set = " << AlignmentOffSet << " Average = " << averageOffset/probeHitCount << endl; */

	}
	std::sort(&LibraryList[0], &LibraryList[Probe.sequenceLength-Probe.NucleotiesInHash-1]);

	int commonValue = 0; 
	int commonValueCount = 0;
	int maxCommonValueCount =0;
	int maxCommonValue =0;
	int zeroCounter = 0; 

	for (int i = 0; i < Probe.sequenceLength - Probe.NucleotiesInHash; i++)
	{
		int value = LibraryList[i] ;
		if (value == 0)
		{
			zeroCounter++; 
			continue; 
		}
		if (commonValue == value)
		{
			commonValueCount++;
			if (commonValueCount > maxCommonValueCount)
			{
				maxCommonValueCount = commonValueCount;
				maxCommonValue = value; 
			}
		} 
		else
		{
			if (commonValueCount > 2)
			{
				//cout << "Found " << commonValue << " a total of " << commonValueCount <<  " times\n";
			}
			commonValue = value; 
			commonValueCount = 0; 
		}
	}

	//maxCommonValue is the result to return
	free (LibraryList);
	timer.Stop (); 
	secondsToAlign = timer.Elapsed (); 
	cout << "total library misses =" << zeroCounter << endl;

	return maxCommonValue;
}//end of 
